In [1]:
from jedi import jedi
from jedi.utils import plot, seedutil, func_generator, init_tools

from functools import partial
import random
import types
import matplotlib.pylab as plt
from matplotlib import animation
from mpl_toolkits.mplot3d import Axes3D
import numpy as np
from __future__ import division
from scipy.integrate import odeint, ode
from numpy import zeros,ones,eye,tanh,dot,outer,sqrt,linspace, \
    cos,pi,hstack,zeros_like,abs,repeat
from numpy.random import uniform,normal,choice
import numpy.linalg as la

from ipywidgets import interact, fixed
from sklearn.decomposition import PCA

%config InlineBackend.figure_format = 'retina'
%matplotlib inline

In [2]:
reload(jedi)
reload(seedutil)


Out[2]:
<module 'jedi.utils.seedutil' from '/Users/simonhaxby/Code/Python/jedi/utils/seedutil.pyc'>

In [3]:
# Setting Seeds
seeds = random.sample(range(100000), 1)

Comparing thresholds


In [4]:
t = linspace(-2,2,1000)
plt.subplot(2,1,1)
plt.plot(t, jedi.step_decode(t));
plt.ylim(-0.5,1.5)

delta = 100
plt.subplot(2,1,2)
plt.plot(t, jedi.sigmoid(delta, t));
plt.ylim(-.5,1.5);


Test Signals

1) Sin Wave


In [5]:
reload(jedi)
reload(seedutil)


Out[5]:
<module 'jedi.utils.seedutil' from '/Users/simonhaxby/Code/Python/jedi/utils/seedutil.pyc'>

In [6]:
# sine-wave target
target = lambda t0: cos(2 * pi * t0/.5)

In [7]:
# simulation parameters for FORCE
dt = .01      # time step
tmax = 10  # simulation length
tstart = 0 # learning start time
tstop = 7  # learning stop time
rho = 1.5   # spectral radius of J
N = 300      # size of stochastic pool
lr = 1.0   # learning rate
pE = .8 # percent excitatory
sparsity = (.1,1,1) # sparsity

In [8]:
noise_mat = np.array([np.random.normal(0,0.3,N) for i in range(1002)])

In [9]:
errors = []
xs = []

for seedling in seeds:
    J, Wz, _, x0, u, w = init_tools.set_simulation_parameters(seedling, N, 1, pE=pE, p=sparsity, rho=rho)
    
    # inp & z are dummy variables
    def model(t0, x, params):
        i = params['index']
        tanh_x = params['tanh_x']
        z = params['z']
        noise = params['noise'][i]
        return (-x + dot(J, tanh_x) + Wz*z + noise)/dt 

    
    x, t, z, _, wu,_ = jedi.force(target, model, lr, dt, tmax, tstart, tstop, x0, w, noise=noise_mat)
    
    xs.append(x)
    error = np.abs(z-target(t))
    errors.append(error)
    
errors = np.array(errors)


Simulation run-time (wall): 5.741 seconds

In [10]:
T = 300
plt.figure(figsize=(12,3))
ax = plt.subplot(111)
ax.spines['right'].set_visible(False)
ax.spines['top'].set_visible(False)
for i in range(20):
    ax.plot(t[:T], x[:T, i], "#696969")
#plt.xlabel('time', fontweight='bold', fontsize=16)



In [11]:
plt.figure(figsize=(12,3))
plot.signal_error(errors, t, tstart, tstop, title= "FORCE (Sin Wave)", burn_in=5)



In [12]:
plt.figure(figsize=(12,4))
plot.target_vs_output_plus_error(t, z, wu, target, offset=1, log=False)



In [13]:
pca = PCA(n_components=2)
pca.fit(x)
pca_x = pca.transform(x).T

In [14]:
reload(plot)


Out[14]:
<module 'jedi.utils.plot' from '/Users/simonhaxby/Code/Python/jedi/utils/plot.pyc'>

In [15]:
tmax = t[-1]
tmin = t[0]
tv = t[:-1]
interact(plot.visualize_2dim_state, time=(tmin, tmax, .1), pca_x=fixed(pca_x), tv=fixed(tv));2


Out[15]:
2

In [16]:
np.save("../data/random/pca_x.npy", pca_x)

In [17]:
pca = PCA(n_components=20)
X_pca = pca.fit_transform(x)

plt.figure(figsize=(8,5))
plt.plot(100*pca.explained_variance_/sum(pca.explained_variance_), 'bo');
plt.xlabel("PCs");
plt.ylabel("Percent of Variance Explained");



In [18]:
tmax = t[-1]
tmin = t[0]
t = t[:-1]

def view_rate(time, xs, t, sim):
    ti = np.argmax(tv >= time)
    plt.title("t=%.2fms" % ti)
    plt.plot(xs[sim][ti, :], np.tanh(xs[sim][ti, :]), 'ro',
             alpha=.1)

In [19]:
interact(view_rate, time=(tmin, tmax, .1), xs=fixed(xs),
         t=fixed(t), sim=fixed(0));



In [20]:
np.save("../data/random/force_sin_xs.npy", np.array(xs))

In [21]:
reload(jedi)
reload(seedutil)


Out[21]:
<module 'jedi.utils.seedutil' from '/Users/simonhaxby/Code/Python/jedi/utils/seedutil.pyc'>

In [22]:
jedi.dforce?

In [23]:
# simulation parameters for FORCE
dt = .01      # time step
tmax = 10  # simulation length
tstart = 0
tstop = 5  # learning stop time
rho = 1.5   # spectral radius of J
N = 300      # size of stochastic pool
lr = 1.0   # learning rate
pE = .8 # percent excitatory
sparsity = (.1,1,1) # sparsity

In [28]:
derrors = []

for seed in seeds:
    J, Wz, _, x0, u, w = init_tools.set_simulation_parameters(seedling, N, 1, pE=pE, p=sparsity, rho=rho)
    
    def model(t0, x, params):
        index = params['index']
        tanh_x = params['tanh_x']
        z = params['z'] 
        noise = params['noise'][index]
        return (-x + dot(J, tanh_x) + Wz*z + noise)/dt 

        x,t,z,w_learn,wu,_ = jedi.dforce(jedi.step_decode, targets, model, lr, dt, tmax, tstart, tstop, x0, w, 
                                     inputs=inputs, noise=noise_mat, pE=pE)


    derror = np.abs(z-target(t))
    derrors.append(derror)
    
derrors = np.array(derrors)


---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-28-b4ba46fc75fe> in <module>()
     15 
     16 
---> 17     derror = np.abs(z-target(t))
     18     derrors.append(derror)
     19 

ValueError: operands could not be broadcast together with shapes (1002,) (1001,) 

In [27]:
plt.figure(figsize=(12,3))
ax = plt.subplot(111)
ax.spines['right'].set_visible(False)
ax.spines['top'].set_visible(False)
for i in range(10):
    ax.plot(t[:T], x[:T, i], "#696969")



In [25]:
plt.figure(figsize=(12,5))
plot.signal_error(derrors, t, tstart, tstop, title="DFORCE (Sin Wave)", burn_in=5)
plot.target_vs_output_plus_error(t, z, wu, target, log=False)


/Users/simonhaxby/anaconda/envs/py27/lib/python2.7/site-packages/numpy/core/_methods.py:59: RuntimeWarning: Mean of empty slice.
  warnings.warn("Mean of empty slice.", RuntimeWarning)
/Users/simonhaxby/anaconda/envs/py27/lib/python2.7/site-packages/numpy/core/_methods.py:70: RuntimeWarning: invalid value encountered in double_scalars
  ret = ret.dtype.type(ret / rcount)
---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
<ipython-input-25-d7a92704a851> in <module>()
      1 plt.figure(figsize=(12,5))
----> 2 plot.signal_error(derrors, t, tstart, tstop, title="DFORCE (Sin Wave)", burn_in=5)
      3 plot.target_vs_output_plus_error(t, z, wu, target, log=False)

/Users/simonhaxby/Code/Python/jedi/utils/plot.pyc in signal_error(errs, t, tstart, tstop, title, burn_in, mean)
     19     if mean:
     20         errs = np.mean(errs, axis=0)
---> 21     ymax = 2*np.max(errs[burn_in:])
     22     plt.plot(t[burn_in:], errs[burn_in:], label="Signal/Output Error")
     23     plt.vlines(tstop,0, ymax, label="Training Stop")

IndexError: invalid index to scalar variable.
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
/Users/simonhaxby/anaconda/envs/py27/lib/python2.7/site-packages/IPython/core/formatters.pyc in __call__(self, obj)
    335                 pass
    336             else:
--> 337                 return printer(obj)
    338             # Finally look for special method names
    339             method = _safe_get_formatter_method(obj, self.print_method)

/Users/simonhaxby/anaconda/envs/py27/lib/python2.7/site-packages/IPython/core/pylabtools.pyc in <lambda>(fig)
    207         png_formatter.for_type(Figure, lambda fig: print_figure(fig, 'png', **kwargs))
    208     if 'retina' in formats or 'png2x' in formats:
--> 209         png_formatter.for_type(Figure, lambda fig: retina_figure(fig, **kwargs))
    210     if 'jpg' in formats or 'jpeg' in formats:
    211         jpg_formatter.for_type(Figure, lambda fig: print_figure(fig, 'jpg', **kwargs))

/Users/simonhaxby/anaconda/envs/py27/lib/python2.7/site-packages/IPython/core/pylabtools.pyc in retina_figure(fig, **kwargs)
    124     """format a figure as a pixel-doubled (retina) PNG"""
    125     pngdata = print_figure(fig, fmt='retina', **kwargs)
--> 126     w, h = _pngxy(pngdata)
    127     metadata = dict(width=w//2, height=h//2)
    128     return pngdata, metadata

/Users/simonhaxby/anaconda/envs/py27/lib/python2.7/site-packages/IPython/core/display.pyc in _pngxy(data)
    607 def _pngxy(data):
    608     """read the (width, height) from a PNG header"""
--> 609     ihdr = data.index(b'IHDR')
    610     # next 8 bytes are width/height
    611     w4h4 = data[ihdr+4:ihdr+12]

AttributeError: 'NoneType' object has no attribute 'index'
<matplotlib.figure.Figure at 0x10cad1f90>

In [ ]:
pca = PCA(n_components=2)
pca.fit(x)

In [ ]:
pca_x = pca.transform(x).T

In [ ]:
tmax = t[-1]
tmin = t[0]
tv = t[:-1]

interact(plot.visualize_2dim_state, time=(tmin, tmax, .1), pca_x=fixed(pca_x), tv=fixed(tv));

In [ ]:
plt.figure(figsize=(7,3))
interact(view_rate, time=(tmin, tmax, .1), xs=fixed(xs),
         t=fixed(t), sim=fixed(0));

In [ ]:
edi.sigmoid(1, np.array(xs))

In [ ]:
pca = PCA(n_components=20)
X_pca = pca.fit_transform(x)

plt.figure(figsize=(8,5))
plt.plot(100*pca.explained_variance_/sum(pca.explained_variance_), 'bo');
plt.xlabel("PCs");
plt.ylabel("Percent of Variance Explained");
plt.title("SFORCE (Sin Wave)");

In [ ]:


In [ ]:


In [ ]:
plt.figure(figsize=(12,3))
plot.cross_signal_error(errors, derrors, t, tstart, tstop, 
                        title="FORCE vs SFORCE (Sin Wave))", burn_in=100)

2) 1D Flip-Flop


In [ ]:
reload(func_generator)

In [ ]:
targets = np.load("../data/stability/flipflop/targets_tmax10sec.npy")
inputs = np.load("../data/stability/flipflop/inputs_tmax10sec.npy")

In [ ]:
# simulation parameters for FORCE
dt = .01      # time step
tmax = 10  # simulation length
tstart = 0
tstop = 5  # learning stop time
rho = 1.02  # spectral radius of J
N = 300      # size of stochastic pool
lr = 1.0   # learning rate
pE = .8 # percent excitatory
sparsity = (.1,1,1) # sparsity
I = 1 # input size

In [ ]:
noise_mat = np.array([np.random.normal(0,.3, N) for i in range(1002)])

In [ ]:
internal_noise = True

errors = []
wus = []
zs = []

for seedling in seeds[:1]:
    J, Wz, Wi, x0, u, w = init_tools.set_simulation_parameters(seedling, N, I, pE=pE, p=sparsity, rho=rho)
    
    # inp & z are dummy variables
    def model(t0, x, params):
        index = params['index']
        tanh_x = params['tanh_x']
        z = params['z']
        inp = params['inputs'][index]
        if internal_noise:
            noise = params['noise'][index]
        else:
            noise = 0
        return (-x + dot(J, tanh_x) + dot(Wi, inp) + Wz*z + noise)/dt
    
    if internal_noise:
        x,t,z,w_learn,wu,_ = jedi.force(targets, model, lr, dt, tmax, tstart, tstop, x0, w, 
                                     inputs=inputs, noise=noise_mat)
    else:
        inputs += noise_mat[:, 0]
        x,t,z,w_learn,wu,_ = jedi.force(targets, model, lr, dt, tmax, tstart, tstop, x0, w, 
                                     inputs=inputs)      

    zs.append(z)
    wus.append(wu)
    
    error = np.abs(z-np.array(targets))
    errors.append(error)
    
errors = np.array(errors)

In [ ]:


In [ ]:
int(tmax/dt+2)

In [ ]:
plt.figure(figsize=(12,5))
plot.signal_error(errors, t, tstart, tstop, title= "FORCE (Flip-Flop)", burn_in=5)

In [ ]:
ind = 0

In [ ]:
plt.figure(figsize=(12,5))
if ind < len(seeds):
    print("Seed: %d" % ind)
    plot.target_vs_output_plus_error(t, zs[ind], wus[ind], targets, offset=0, log=False)
    ind+=1

In [ ]:
T = 500
plt.figure(figsize=(12,3))
ax = plt.subplot(111)
ax.spines['right'].set_visible(False)
ax.spines['top'].set_visible(False)
for i in range(10):
    ax.plot(t[:T], x[:T, i], "#696969")
#plt.xlabel('time', fontweight='bold', fontsize=16)

In [ ]:
plt.figure(figsize=(12,5))
plt.legend()
plt.subplot(2,1,2)
for i in range(10):
    plt.plot(t[:], x[:, i]);
plt.subplot(2,1,1)
plt.plot(t, np.array(z), lw=2, label="output");

In [ ]:
pca = PCA(n_components=20)
X_pca = pca.fit_transform(x)

plt.figure(figsize=(8,5))
plt.plot(100*pca.explained_variance_/sum(pca.explained_variance_), 'bo');
plt.xlabel("PCs");
plt.ylabel("Percent of Variance Explained");
plt.title("FORCE (Flip Flop)");

In [ ]:
pca = PCA(n_components=2)
pca.fit(x)
pca_x = pca.transform(x).T

In [ ]:
tmax = t[-1]
tmin = t[0]
tv = t[:-1]

interact(plot.visualize_2dim_state, time=(tmin, tmax, .1), pca_x=fixed(pca_x), tv=fixed(tv));

In [ ]:
np.save("../data/random/pca_x_ff_noise.npy", pca_x)

In [ ]:
"({0}, {1})".format(.1,.2)

In [ ]:
pca_x.shape

In [ ]:
# simulation parameters for FORCE
dt = .01      # time step
tmax = 10  # simulation length
tstop = 5  # learning stop time
rho = 1.02   # spectral radius of J
N = 300      # size of stochastic pool
lr = 1.0   # learning rate
pE = .8 # percent excitatory
sparsity = (.1,1,1) # sparsity
I = 1 # input size

In [ ]:
internal_noise = True

derrors = []
wus = []
zs = []

for seedling in seeds[:2]:
    J, Wz, Wi, x0, u, w = init_tools.set_simulation_parameters(seedling, N, I, pE=pE, p=sparsity, rho=rho)
    
    # inp & z are dummy variables
    def model(t0, x, params):
        index = params['index']
        tanh_x = params['tanh_x']
        z = params['z']
        inp = params['inputs'][index]
        if internal_noise:
            noise = params['noise'][index]
        else:
            noise = 0
        return (-x + dot(J, tanh_x) + dot(Wi, inp) + Wz*z + noise)/dt

    if internal_noise:
        x,t,z,w_learn,wu,_ = jedi.dforce(jedi.step_decode, targets, model, lr, dt, tmax, tstart, tstop, x0, w, 
                                     inputs=inputs, noise=noise_mat, pE=pE)
    else:
        inputs += noise_mat[:, 0]
        x,t,z,w_learn,wu,_ = jedi.dforce(jedi.step_decode, targets, model, lr, dt, tmax, tstart, tstop, x0, w, 
                                     inputs=inputs, pE=pE)      

    zs.append(z)
    wus.append(wu)
    
    derror = np.abs(z-np.array(targets))
    derrors.append(derror)
    
derrors = np.array(derrors)

In [ ]:
T = 500
plt.figure(figsize=(12,3))
ax = plt.subplot(111)
ax.spines['right'].set_visible(False)
ax.spines['top'].set_visible(False)
for i in range(10):
    ax.plot(t[:T], x[:T, i], '#696969')
#plt.xlabel('time', fontweight='bold', fontsize=16)

In [ ]:
plt.figure(figsize=(12,5))
plot.signal_error(derrors, t, tstart, tstop, title="SFORCE (Flip-Flop)", burn_in=5)

In [ ]:
# Setting seed index
ind = 0

In [ ]:
plt.figure(figsize=(12,5))
if internal_noise:
    plot.target_vs_output_plus_error(t, zs[ind], wus[ind], targets, offset=0, log=False)
else:
    plot.target_vs_output_plus_error(t, zs[ind], wus[ind], targets+noise_mat[:,0], offset=0, log=False)

In [ ]:
plt.figure(figsize=(12,4))
plt.legend()
plt.subplot(2,1,2)
for i in range(20):
    plt.plot(t[:], x[:, i]);
plt.subplot(2,1,1)
plt.plot(t, np.array(z), lw=2, label="output");

In [ ]:
pca = PCA(n_components=20)
X_pca = pca.fit_transform(x)

plt.figure(figsize=(8,5))
plt.plot(100*pca.explained_variance_/sum(pca.explained_variance_), 'bo');
plt.xlabel("PCs");
plt.ylabel("Percent of Variance Explained");
plt.title("DFORCE (Flip Flop)");

In [ ]:
pca = PCA(n_components=3)
pca.fit(x)
pca_x = pca.transform(x).T

In [ ]:
tmax = t[-1]
tmin = t[0]
tv = t[:-1]

interact(plot.visualize_2dim_state, time=(tmin, tmax, .1), pca_x=fixed(pca_x), tv=fixed(tv));

In [ ]:
interact(plot.visualize_3dim_state, time=(tmin, tmax, .1), pca_x=fixed(pca_x), tv=fixed(tv));

In [ ]:
import seaborn

In [ ]:


In [ ]:
# cross mean signal error
plt.figure(figsize=(12,5))
plot.cross_signal_error(errors, derrors, t, tstart, tstop, 
                        title="FORCE vs SFORCE (Flip-Flop)", burn_in=5)

In [ ]:
# Setting seed index
ind = 0

In [ ]:
errors[0]

In [ ]:
plt.figure(figsize=(12,5))
print("Seed: %d" % ind)
plot.cross_signal_error(errors[ind], derrors[ind], t, tstart, tstop, 
                        title="FORCE vs SFORCE, (Flip-Flop)", 
                        burn_in=5, mean=False)
ind+=1

3) Lorentz Attractor


In [ ]:
# Parameters specified by Abbott 2009.
def lorentz((x, y, z), t0, sigma=10., beta=8./3, rho=28.0):
    return [sigma * (y - x), x * (rho - z) - y, x * y - beta * z]

In [ ]:
break_in = 500
T = 1501 # period
x0 = np.random.randn(3) # starting vector
t_= np.linspace(0, 60, T)
lorenz = odeint(lorentz, x0, t_)/20
targets = lorenz[break_in:,0]
t_in = t_[break_in:]

In [ ]:
# Visualizing Lorentz attractor
plt.figure(figsize=(12,3))
plt.plot(t_in, targets);
plt.xlim([min(t_in), max(t_in)])
plt.ylim([-1.5,1.5])
plt.xlabel('time', fontsize=14);
plt.ylabel('y', fontsize=14);

In [ ]:
# simulation parameters for FORCE
dt = .01      # time step
tmax = 10  # simulation length
tstop = 8  # learning stop time
rho = 1.2   # spectral radius of J
N = 1000     # size of stochastic pool
lr = 1.0   # learning rate
pE = None # percent excitatory
sparsity = (.1,1,1) # sparsity

In [ ]:
errors = []
wus = []
zs = []

for seedling in seeds:
    J, Wz, _, x0, u, w = init_tools.set_simulation_parameters(seedling, N, 1, pE=pE, p=sparsity, rho=rho)
    
    # inp & z are dummy variables
    def model(t0, x, tanh_x, inp, z): 
        return (-x + dot(J, tanh_x) + Wz*z)/dt 
    
    x, t, z, _, wu,_ = jedi.force(targets, model, lr, dt, tmax, tstart, tstop, x0, w, 0)

    zs.append(z)
    wus.append(wu)
    
    error = np.abs(z[1:]-np.array(targets))
    errors.append(error)
    
errors = np.array(errors)

In [ ]:
plt.figure(figsize=(12,7))
ind = 0  
tstart = 1
plot.target_vs_output_plus_error(t[tstart:], zs[ind][tstart:], 
                                 wus[ind][tstart:], targets[tstart:], offset=1, ylim=[-1,2])

In [ ]:
# Visualizing activities of first 20 neurons
T = 500
K = 20

plt.figure(figsize=(12,4))
plt.subplot(211)
plt.title("Neuron Dynamics");
for i in range(K):
    plt.plot(t[:T], x[:T, i]);
    
plt.subplot(212)
for i in range(K):
    plt.plot(t[-T:], x[-T:, i]);
    plt.xlim(t[-T], t[-1])

In [ ]:
plt.figure(figsize=(12,5))
plot.signal_error(errors, t[1:], tstop, title= "FORCE", burn_in=5)

In [ ]:
ind = 0

In [ ]:
pca = PCA(n_components=20)
X_pca = pca.fit_transform(x)

plt.figure(figsize=(8,5))
plt.plot(100*pca.explained_variance_/sum(pca.explained_variance_), 'bo');
plt.xlabel("PCs");
plt.ylabel("Percent of Variance Explained");
plt.title("FORCE (Flip Flop)");

In [ ]:
pca = PCA(n_components=3)
pca.fit(x)
pca_x = pca.transform(x).T

In [ ]:
tmax = t[-1]
tmin = t[0]
tv = t[:-1]

interact(plot.visualize_3dim_state, time=(tmin, tmax, .1), pca_x=fixed(pca_x), tv=fixed(tv));

In [ ]:


In [ ]:
# simulation parameters for FORCE
dt = .01      # time step
tmax = 10  # simulation length
tstop = 8  # learning stop time
rho = 1.02   # spectral radius of J
N = 1000     # size of stochastic pool
lr = 1.0   # learning rate
pE = None # percent excitatory
sparsity = (.1,1,1) # sparsity

In [ ]:
derrors = []
wus = []
zs = []

for seedling in seeds:
    J, Wz, _, x0, u, w = init_tools.set_simulation_parameters(seedling, N, 1, pE=pE, p=sparsity, rho=rho)
    
    # inp & z are dummy variables
    def model(t0, x, tanh_x, inp, z): 
        return (-x + dot(J, tanh_x) + Wz*z)/dt 
    
    x, t, z, _, wu,_ = jedi.dforce(jedi.step_decode, targets, model, lr, dt, tmax, tstart, tstop, x0, w, 0, pE=.8)

    zs.append(z)
    wus.append(wu)
    
    derror = np.abs(z[1:]-np.array(targets))
    derrors.append(derror)
    
derrors = np.array(derrors)

In [ ]:
plt.figure(figsize=(12,5))
ind = 0 
plot.target_vs_output_plus_error(t, zs[ind], wus[ind], targets, offset=1, log=False)

In [ ]:
plt.figure(figsize=(12,4))
plot.signal_error(derrors, t[1:], tstop, title= "SFORCE", burn_in=5)

In [ ]:
plt.figure(figsize=(12,5))
print("Seed: %d" % ind)
plot.signal_error(derrors[ind], t[1:], tstop, 
                  title= "SFORCE", burn_in=5, mean=False)
plot.target_vs_output_plus_error(t, zs[ind], wus[ind], targets, offset=1, log=False)
ind+=1

In [ ]:
pca = PCA(n_components=20)
X_pca = pca.fit_transform(x)

plt.figure(figsize=(8,5))
plt.plot(100*pca.explained_variance_/sum(pca.explained_variance_), 'bo');
plt.xlabel("PCs");
plt.ylabel("Percent of Variance Explained");
plt.title("SFORCE (Flip Flop)");

In [ ]:
pca = PCA(n_components=3)
pca.fit(x)
pca_x = pca.transform(x).T

In [ ]:
tmax = t[-1]
tmin = t[0]
tv = t[:-1]

interact(plot.visualize_3dim_state, time=(tmin, tmax, .1), pca_x=fixed(pca_x), tv=fixed(tv));

In [ ]:


In [ ]:
# cross mean signal error
plot.cross_signal_error(errors, derrors, t[1:], tstop, 
                        title="FORCE vs SFORCE (Lorentz)", burn_in=5)

In [ ]: